{
  "cells": [
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "%matplotlib inline"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "\n# Curve Fitting with Bayesian Ridge Regression\n\n\nComputes a Bayesian Ridge Regression of Sinusoids.\n\nSee `bayesian_ridge_regression` for more information on the regressor.\n\nIn general, when fitting a curve with a polynomial by Bayesian ridge\nregression, the selection of initial values of\nthe regularization parameters (alpha, lambda) may be important.\nThis is because the regularization parameters are determined by an iterative\nprocedure that depends on initial values.\n\nIn this example, the sinusoid is approximated by a polynomial using different\npairs of initial values.\n\nWhen starting from the default values (alpha_init = 1.90, lambda_init = 1.),\nthe bias of the resulting curve is large, and the variance is small.\nSo, lambda_init should be relatively small (1.e-3) so as to reduce the bias.\n\nAlso, by evaluating log marginal likelihood (L) of\nthese models, we can determine which one is better.\nIt can be concluded that the model with larger L is more likely.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "print(__doc__)\n\n# Author: Yoshihiro Uchida <nimbus1after2a1sun7shower@gmail.com>\n\nimport numpy as np\nimport matplotlib.pyplot as plt\n\nfrom sklearn.linear_model import BayesianRidge\n\n\ndef func(x): return np.sin(2*np.pi*x)\n\n\n# #############################################################################\n# Generate sinusoidal data with noise\nsize = 25\nrng = np.random.RandomState(1234)\nx_train = rng.uniform(0., 1., size)\ny_train = func(x_train) + rng.normal(scale=0.1, size=size)\nx_test = np.linspace(0., 1., 100)\n\n\n# #############################################################################\n# Fit by cubic polynomial\nn_order = 3\nX_train = np.vander(x_train, n_order + 1, increasing=True)\nX_test = np.vander(x_test, n_order + 1, increasing=True)\n\n# #############################################################################\n# Plot the true and predicted curves with log marginal likelihood (L)\nreg = BayesianRidge(tol=1e-6, fit_intercept=False, compute_score=True)\nfig, axes = plt.subplots(1, 2, figsize=(8, 4))\nfor i, ax in enumerate(axes):\n    # Bayesian ridge regression with different initial value pairs\n    if i == 0:\n        init = [1 / np.var(y_train), 1.]  # Default values\n    elif i == 1:\n        init = [1., 1e-3]\n        reg.set_params(alpha_init=init[0], lambda_init=init[1])\n    reg.fit(X_train, y_train)\n    ymean, ystd = reg.predict(X_test, return_std=True)\n\n    ax.plot(x_test, func(x_test), color=\"blue\", label=\"sin($2\\\\pi x$)\")\n    ax.scatter(x_train, y_train, s=50, alpha=0.5, label=\"observation\")\n    ax.plot(x_test, ymean, color=\"red\", label=\"predict mean\")\n    ax.fill_between(x_test, ymean-ystd, ymean+ystd,\n                    color=\"pink\", alpha=0.5, label=\"predict std\")\n    ax.set_ylim(-1.3, 1.3)\n    ax.legend()\n    title = \"$\\\\alpha$_init$={:.2f},\\\\ \\\\lambda$_init$={}$\".format(\n            init[0], init[1])\n    if i == 0:\n        title += \" (Default)\"\n    ax.set_title(title, fontsize=12)\n    text = \"$\\\\alpha={:.1f}$\\n$\\\\lambda={:.3f}$\\n$L={:.1f}$\".format(\n           reg.alpha_, reg.lambda_, reg.scores_[-1])\n    ax.text(0.05, -1.0, text, fontsize=12)\n\nplt.tight_layout()\nplt.show()"
      ]
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.6.9"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}